In [1]:
%pylab inline
import numpy as np # convert list to array
import matplotlib.pyplot as plt
from ipywidgets import *
#from scipy import special
#import scipy.special as sp
from scipy.special import jn,fresnel
%pylab is deprecated, use %matplotlib inline and import the required libraries. Populating the interactive namespace from numpy and matplotlib
Talbot carpet¶
The wave from a phase grating producing a sinusoidal wavefront with amplitude H is:
$ \psi(x,z) =\sum_{n=-\infty}^\infty i^n J_n(\frac{2H}{\lambda}) \, \exp\left[ i\pi \left(2n\frac{x}{a}-n^2 \frac{z}{z_T}\right)\right] = J_0(\frac{2H}{\lambda}) + 2\sum_{n=1}^\infty i^n J_n(\frac{2H}{\lambda}) \, \cos(2 \pi n\frac{x}{a}) \exp\left( -i\pi n^2 \frac{z}{z_T}\right), $
where the $J_n$ are Bessel functions.
Talbot period:¶
$\psi(x,z) = \psi(x+\frac{a}{2},z+z_T)$
Therefore the intensity is indeed periodic with the Talbot period, apart from a shift of $a/2$ between successive images.
Papers:¶
Berry, M V, Marzoli, I and Schleich, W, 2001 ‘Quantum carpets, carpets of light’ Physics World (June), 39-44.¶
and
Berry, M V and Bodenschatz, E, 1999 ‘Caustics, multiply-reconstructed by Talbot interference’, J.Mod.Opt 46 349 – 365.¶
In [2]:
def psi(x,z,param):
# Nmax --> n=1 to n=Nmax in the sum
H,lam,Nmax = param
### NEM JO, ha nvec egy array es x, z -t meshgrid-del hivjuk az nvec*x, es nvec*z miatt. lista*matrix =???
#nvec=array(range(1, Nmax+1))
#tmp=(1j)**nvec*jn(nvec,2*H/lam)*cos(2*pi*nvec*x/a)*exp(-1j*pi*nvec**2*z/zT)
#res=jn(0,2*H/lam)+2*sum(tmp)
## Ezert a for loop a jo megoldas!!!!
sum=0
for n in range(1, Nmax+1):
tmp= (1j)**n*jn(n,2*pi*H/lam)*cos(2*pi*n*x)*exp(-1j*pi*n**2*z)
sum+=tmp
sum= jn(0,2*pi*H/lam)+2*sum
return(sum)
In [3]:
def psi_orig(x,z,param):
# for loop, de az osszegzesben: Nmax --> n=-Nmax to n=Nmax
# NEM ezt hasznaljuk!!!
H,lam,Nmax = param
sum=0
for n in range(-Nmax, Nmax+1):
tmp= (1j)**n*jn(n,2*pi*H/lam)*exp(1j*pi*(2*n*x-n**2*z))
sum+=tmp
return(sum)
In [4]:
'''
rajzolashoz, $\psi$ hivja.
'''
def wave_field(xr,yr,param,Np=100):
#mintavételezési pontok legyártása
x2,y2 = meshgrid(linspace(-xr,xr,Np),linspace(0,yr,Np))
H,lam,Nmax = param
z2u=psi(x2,y2,param)
z2=abs(z2u)**2
return(x2,y2,z2)
In [5]:
param=[400*pi,1,17]
psi_orig(0.27,0.57,param), psi(0.27,0.57,param)
Out[5]:
((0.0169148410475245-0.021710095379094452j), (0.016914841047523957-0.021710095379094442j))
In [6]:
xlim=0.75 ## plotted range in x
zlim=2.25 ## plotted range in y (a fuggveny definicojaban)
Np=150 # a mintavetelezesi pontok, Np x Np szamu racspontban szamoljuk a psi-t
#xlim=0.1 ## plotted range in x
#zlim=0.025 ## plotted range in y (a fuggveny definicojaban)
#Np=200 # a mintavetelezesi pontok, Np x Np szamu racspontban szamoljuk a psi-t
### parameters:
a=1 ## racsallando
#H=400*pi
H= 1/(1/50*4*pi**2) ## a szinusz amplitudoja
#H= 40/pi**2 ## a szinusz amplitudoja
lam =1 ## hullamhossz
Nmax = int(2*2*H/lam) ##az osszegzesben a tagok szama
print('2*H/lam =',round(2*H/lam,2),', N_max = ', Nmax)
param= [H,lam,Nmax]
zT=a**2/lam ## Talbot-periodus z-ben
rho=1/(4*pi**2)*lam/H
print('z_T, Talbot period along z = ' ,zT, "\n")
print('rho = ' ,rho)
print('dx/a = rho = ' ,rho)
print('dz/z_T = ' ,8*pi*rho**2/sqrt(1+(4*pi*rho)**2),"\n")
print('cusp: dx/a = ' ,rho**(3/4))
print('cusp: dz/a = ' ,rho**(3/2),"\n")
x2,z2,fv=wave_field(xlim,zlim,param,Np)
fvmax=max(flatten(fv))
print('|psi_max|^2 =',fvmax)
### plotting
figsize(10,8)
#figsize(xfig_meret,yfig_meret)
res=contourf(x2,z2,fv,levels=linspace(0.,fvmax,100),cmap='rainbow') # cmap='viridis', cmap='Reds' cmap='gray'
#axes().set_aspect('equal')
#res=contourf(x2,y2,z2,cmap='rainbow')
title('Talbot carpet',fontsize=16);
xlabel(r'$x/a$',fontsize=16)
ylabel(r'$z/z_T$',fontsize=16);
#colorbar();
colorbar(res, fraction=0.046, pad=0.04);
2*H/lam = 2.53 , N_max = 5 z_T, Talbot period along z = 1.0 rho = 0.02 dx/a = rho = 0.02 dz/z_T = 0.009749883354388342 cusp: dx/a = 0.053182958969449884 cusp: dz/a = 0.00282842712474619 |psi_max|^2 = 2.7520515067080926
In [7]:
xlim=0.4 ## plotted range in x
zlim=0.07 ## plotted range in y (a fuggveny definicojaban)
Np=200 # a mintavetelezesi pontok, Np x Np szamu racspontban szamoljuk a psi-t
#xlim=0.1 ## plotted range in x
#zlim=0.025 ## plotted range in y (a fuggveny definicojaban)
#Np=200 # a mintavetelezesi pontok, Np x Np szamu racspontban szamoljuk a psi-t
### parameters:
a=1 ## racsallando
#H=400*pi
H= 1/(1/100*4*pi**2) ## a szinusz amplitudoja
#H= 40/pi**2 ## a szinusz amplitudoja
lam =1 ## hullamhossz
Nmax = int(2*2*H/lam) ##az osszegzesben a tagok szama
print('2*H/lam =',round(2*H/lam,2),', N_max = ', Nmax)
param= [H,lam,Nmax]
zT=a**2/lam ## Talbot-periodus z-ben
rho=1/(4*pi**2)*lam/H
print('z_T, Talbot period along z = ' ,zT, "\n")
print('rho = ' ,rho)
print('dx/a = rho = ' ,rho)
print('dz/z_T = ' ,8*pi*rho**2/sqrt(1+(4*pi*rho)**2),"\n")
print('cusp: dx/a = ' ,rho**(3/4))
print('cusp: dz/a = ' ,rho**(3/2),"\n")
x2,z2,fv=wave_field(xlim,zlim,param,Np)
fvmax=max(flatten(fv))
print('|psi_max|^2 =',fvmax)
### plotting
figsize(10,8)
#figsize(xfig_meret,yfig_meret)
res=contourf(x2,z2,fv,levels=linspace(0.,fvmax,100),cmap='rainbow') # cmap='viridis', cmap='Reds' cmap='gray'
#axes().set_aspect('equal')
#res=contourf(x2,y2,z2,cmap='rainbow')
title('Talbot carpet',fontsize=16);
xlabel(r'$x/a$',fontsize=16)
ylabel(r'$z/z_T$',fontsize=16);
#colorbar();
colorbar(res, fraction=0.046, pad=0.04);
2*H/lam = 5.07 , N_max = 10 z_T, Talbot period along z = 1.0 rho = 0.01 dx/a = rho = 0.01 dz/z_T = 0.002493662078269732 cusp: dx/a = 0.03162277660168379 cusp: dz/a = 0.001 |psi_max|^2 = 5.409228101385689
In [ ]: